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ABSTRACT 


Efforts  to  update  eurrent  wave  speed  models  of  the  Middle  East  require  a  thoroughly  tested  database  of  sources  and 
recordings.  Recordings  of  seismie  waves  traversing  the  region  from  Tibet  to  the  Red  Sea  will  be  the  principal  metrie 
in  guiding  improvements  to  the  current  wave  speed  model.  Preeise  characterizations  of  the  earthquakes,  specifically 
depths  and  faulting  meehanisms,  are  essential  to  avoid  mapping  souree  errors  into  the  refined  wave  speed  model. 
Errors  associated  with  the  source  are  manifested  in  amplitude  and  phase  changes.  Souree  depths  and  paths  near 
nodal  planes  are  particularly  error  prone  as  small  changes  may  severely  affect  the  resulting  wavefield.  Onee  sourees 
are  quantified,  regions  requiring  refinement  will  be  highlighted  using  adjoint  tomography  methods  based  on  spectral 
element  simulations  (Komatitsch  and  Tromp,  1999). 

An  initial  daiabase  of  250  regional  Middle  Eastern  events  from  1990-2007,  was  inverted  for  depth  and  focal 
mechanism  using  teleseismie  arrivals  (Kikuehi  and  Kanamori,  1982)  and  regional  surfaee  and  body  waves  (Zhao 
and  Helmberger,  1994).  From  this  initial  database,  we  reinterpreted  a  large,  well-reeorded  subset  of  201  events 
through  a  direet  comparison  between  data  and  synthetics  based  upon  a  centroid  moment  tensor  inversion  (Liu  et  al., 
2004).  Evaluation  was  done  using  both  a  ID  reference  model  (Dziewonski  and  Anderson,  1981)  at  periods  greater 
than  80  seconds  and  a  3D  model  (Kustowski  et  al.,  2008)  at  periods  of  25  seeonds  and  longer.  The  final  souree 
reinterpretations  will  be  within  the  3D  model,  as  this  is  the  initial  starting  point  for  the  adjoint  tomography. 
Transitioning  from  a  ID  to  3D  wave  speed  model  shows  dramatic  improvements  when  comparisons  are  done  at 
shorter  periods,  e.g.,  25  s.  Synthetics  from  the  ID  model  were  created  through  mode  summations  while  those  from 
the  3D  simulations  were  created  using  the  spectral  element  method. 

To  further  assess  errors  in  source  depth  and  focal  mechanism,  comparisons  between  the  three  methods  were  made. 
These  comparisons  help  to  identify  problematic  stations  and  sources  whieh  may  bias  the  final  solution  Estimates  of 
standard  errors  were  generated  for  each  event’s  source  depth  and  focal  mechanism  to  identify  poorly  constrained 
events.  A  final,  wcll-charactcrizcd  set  of  sourees  and  stations  will  be  used  to  iteratively  improve  the  wave  speed 
model  of  the  Middle  East.  After  a  few  iterations  during  the  adjoint  inversion  process,  the  sourees  will  be  reexamined 
and  relocated  to  further  reduee  mapping  of  source  errors  into  structural  features. 

Finally,  efforts  continue  in  developing  the  infrastructure  required  to  “quickly”  generate  event  kernels  at  the  /7th 
iteration  and  invert  for  a  new,  (tf+/)th,  wave  speed  model  of  the  Middle  East.  While  development  of  the 
infrastructure  proceeds,  initial  tests  using  a  limited  number  of  events  demonstrate  that  the  3D  model,  while  showing 
vast  improvement  compared  to  the  ID  model,  still  requires  substantial  modifications.  Employing  our  new,  full 
source  set  and  iterating  the  adjoint  inversions  at  successively  shorter  periods  will  lead  to  significant  ehanges  and 
refined  wave  speed  structures  of  the  Middle  East. 
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OBJECTIVES 


Improved  models  of  3D  wave  speeds  within  the  greater  Middle  East, 
from  the  Turkish  Plateau  to  the  eastern  edge  of  Tibet  (Figure  1),  will 
increase  our  ability  to  discriminate  between  natural  and  man-made 
events,  locate  these  events,  identify  source  depths,  and  determine 
magnitudes.  Current  wave  speed  models  of  the  Middle  East  will  be 
improved  through  an  adjoint  tomography  method  to  iteratively  update 
the  models  in  three  dimensions.  The  initial  step  towards  any  well- 
constrained  tomography  model  involves  the  compilation  of  a  database 
of  relocated  and  characterized  set  of  sources  and  waveforms.  Events 
have  been  re-inverted  in  1 D  and  3D  and  agree  fairly  well  with 
currently  available  solutions.  The  adjoint  tomography  method  (Tromp 
et  al.,  2005)  uses  full  seismic  waveforms  as  a  measure  of  misfit  of  the 
current  model  iteration.  Using  the  differences  between  data  and 
synthetics,  adjoint  sources  are  constructed  and  used  to  generate  kernels 
required  to  update  the  current  iteration  model.  In  future  years,  an 
improved  wave  speed  model  of  the  Middle  East  will  be  presented  and 
distributed,  built  on  the  foundation  of  the  seismic  waveform  data  set. 

RESEARCH  ACCOMPLISHED 


Figure  1.  Map  of  the  available  sources 

and  stations  in  the  Middle  East 
for  the  adjoint  tomographic 
inversion. 


Tomographic  updates  will  use  starting  model  S2.9EA  (Kustowski  et  al.,  2008), 
but  to  avoid  mapping  source  errors  into  the  wave  speed  model,  sources  need  to 
be  re-evaluated  within  this  model.  Re-evaluation  includes  an  assessment  of  the 
faulting  parameters  and  locations  using  full  waveforms.  Multiple  evaluations 
were  conducted  using  synthetic  and  real  waveform  data  to  identify  and 
eliminate  possible  problems  and  sources  of  error. 

Synthetic  Event  Characterization 

Using  a  synthetic  waveform  data  set,  the  capabilities  of  the  centroid  moment 
tensor  (CMT)  inversion  (Liu  et  al.,  2004)  were  assessed.  Figure  2  displays 
results  from  different  parameterized  inversions  in  an  attempt  to  recover  the 
solution  in  red  at  the  top.  Synthetics  for  the  optimal  solution  and  subsequent 
inversions  were  computed  in  model  S2.9EA  (Kustowski  et  al.,  2008). 
Synthetics  for  the  initial  starting  point  were  different  in  hypocenter  and 
moment  tensor  elements  to  force  the  CMT  inversion  to  recover  the  optimal 
solution.  A  simple  CMT  inversion,  6  parameters  (black,  top  row),  that  does  not 
solve  for  location  or  depth  does  a  reasonable  job  in  solution  recovery.  As  more 
parameters  are  added,  including  depth,  7  parameters  (black,  middle  row),  and 
epicenter,  9  parameters  (black,  bottom  row),  the  variance  is  further  reduced  as 
the  original  moment  tensor  and  depth  are  well  represented.  The  left  and  right 
columns  of  Figure  2  show  inversions  for  only  double  couple  (DC)  or  zero-trace 
(ZT)  solutions.  Zero-trace  solutions  perform  better  at  representing  the  source 
as  indicated  by  the  smaller  variance  across  the  range  of  parameters  used,  6-9. 
The  new  database  of  events  will  be  built  using  the  zero-trace  (ZT)  solutions. 


ID  Event  Characterization 
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Figure  2.  Performance  of  the  CMT 
inversion  to  recover  faulting 
parameters  and  source  depth 
from  a  synthetic  event  (red) 
using  several  different 
parameterizations  (black). 
Variance  reductions  and  optimal 
depth  for  each  parameterization 
are  displayed  beneath  each 
moment  tensor.  See  text  for 
further  explanation. 


Using  the  CMT  inversion  methodology  of  Liu  et  al.  (2004),  -200  events  in  the  Middle  East  were  re-characterized 
using  the  ID  PREM  model  (Dziewonski  and  Anderson,  1981)  at  periods  of  80  seconds  and  longer  to  accomplish 
two  goals.  First,  the  procedure  of  processing  large  amounts  of  data,  comparing  these  to  synthetics,  and  reevaluating 
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the  source  parameters  and  locations  needed  to  be  assessed  and  streamlined.  This  procedure  needs  to  be 
straightforward  and  able  to  avoid  problematic  areas,  as  it  will  be  used  in  the  initial  3D  event  re-evaluation  as  well  as 
two  to  three  more  times  during  the  adjoint  inversions.  Second,  problematic  events  and  stations  needed  to  be 
identified  and  removed  before  building  the  final  waveform  data  set.  Problematic  stations  and  waveforms  with 
dropouts  and  poorly  characterized  amplitude  responses  can  negatively  influence  a  CMT  inversion.  These  stations 
were  removed  before  any  3D  CMT  inversions  were  performed.  Solutions  from  the  long  period  ID  solutions  are 
shown  in  Figure  3  (left).  Faulting  parameters  and  depths  agree  reasonably  well  with  those  from  The  Global  CMT 
Project.  Depths  from  the  ID  model  trend  towards  shallower  depths,  which  is  most  likely  due  to  the  3D  model 
corrections  used  at  The  Global  CMT  Project. 

1 D  Solutions  3D  Solutions 


Figure  3.  Comparison  of  new  solutions  using  a  ID  model,  left,  and  3D  model,  right.  Fault  plane  solutions  do 
not  change  appreciably  when  using  either  model,  but  the  depths  do  change,  especially  when  the 
events  are  deep.  Catalog  using  the  ID  model  includes  -200  events,  and  the  3D  model  catalog 
currently  has  -100  events  (50%  finished). 


3D  Event  Characterization 

The  same  events  were  also  re-evaluated  within  the  S2.9EA  model  (Kustowski  et  al.,  2008)  at  a  much  shorter  period 
of  25  seconds.  Figure  3  (right)  shows  the  current  set  of  solutions,  -100  events  (50%  complete),  inverted  within  the 
3D  wave  speed  model.  Faulting  solutions  match  well  with  the  published  solutions,  and  in  contrast  to  the  1 D 
solutions,  the  depths  agree  fairly  well  with  published  depths.  While  the  solutions  and  depths  match  reasonably  well 
when  using  either  ID  or  3D  wave  speed  models,  the  3D  wave  speed  synthetics  use  a  much  smaller  minimum  period, 
25  seconds,  that  incorporates  complicated  body  wave  arrivals  and  more  pronounced  surface  wave  arrivals  and 
dispersion. 

Figure  4  displays  vertical  component  data,  black,  and  synthetics,  red,  for  a  single  event  on  2005/03/13  in  the 
southern  Zagros  Mountains  recorded  by  two  stations  to  the  north.  Synthetics  from  the  1 D  PREM  model  at  long  and 
short  periods  are  displayed  in  the  first  two  columns.  At  long  periods,  the  data  and  synthetics  arc  relatively  simple, 
with  identifiable,  small  amplitude  body  wave  arrivals  and  much  larger  amplitude,  fundamental  mode,  surface  wave 
arrivals.  Synthetics  match  the  absolute  arrivals  times  rather  well,  but  the  amplitudes  arc  misfit  for  each  station  at  the 
longest  period,  80  seconds.  As  shorter  periods  are  included,  middle  column  at  25  seconds,  the  synthetics  remain 
simple  arrivals  as  the  data  starts  to  shows  dispersion  due  to  traversing  the  continent,  see  800  to  1  100  seconds  for 
station  KURK.II.LHZ  and  1000  to  1400  seconds  for  station  KIEV.IU.LHZ.  Body  wave  arrivals  also  show  more 
complex  propagation  as  compared  to  the  ID  synthetics,  including  larger  amplitudes,  longer  duration  coda,  600  to 
700  seconds  for  KURK.II.LHZ  and  700  to  900  for  KIEV.IU.LHZ  in  Figure  4  middle  column.  Simulations  using  the 
S2.9EA  model  (Kustowski  et  al.,  2008)  model  reproduce  the  recorded  data  more  closely  than  the  ID  model  at  25 
seconds.  Figure  4  right  column.  Use  of  the  3D  wave  speed  model  provides  a  more  robust,  shorter  period  estimation 
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of  the  faulting  parameters  without  losing  the  essential,  long  period  estimate  of  the  moment/magnitude.  The  more 
complex  propagation  paths,  represented  by  the  waveform  coda  in  the  data  and  now  the  synthetics,  do  not  influence 
the  estimates  of  source  epicenter,  faulting  parameters,  and  more  importantly,  source  depth. 


Figure  4.  Comparison  of  vertical  data  and  synthetics  at  long  periods,  80  seconds,  and  short  periods,  25 

seconds,  using  ID  and  3D  wave  speed  models.  From  the  event  shown  in  the  map  on  the  left,  the  data 
from  two  stations  arc  compared  to  the  different  models.  Inclusion  of  the  3D  wave  speed  model 
improves  the  match  between  the  data  and  synthetics,  especially  for  the  surface  wave  dispersion  and 
the  complex  body  wave  arriv  als. 


Wave  Speed  Model  Performance 

As  displayed  in  Figure  4,  right  column,  the  complexity  of  the  waveforms  is  well  represented  by  the  3D  wave  speed 
model,  including  the  complex  body  wave  arrivals  and  continental  generated  dispersion  of  the  surface  waves.  While 
the  general  character  of  the  wavefield  is  reproduced  by  the  synthetics,  the  small  variations  in  travel  time  and 
amplitude  will  be  used  to  update  the  current  iteration  wave  speed  model  within  the  adjoint  tomography.  Minor 
mismatches  between  data  and  synthetics  are  visible  in  a  comparison  for  a  single  event  on  2005/02/22  recorded  at 
stations  around  the  Middle  East  (Figure  5).  Data  from  the  tangential  component,  black,  is  compared  with  3D 
synthetics,  red,  and  plotted  against  azimuth  to  highlight  different  propagation  paths.  Paths  towards  the  north  show 
large  amplitude,  dispersive  Love  waves  that  are  reasonably  matched  by  the  simulations,  specifically  the  dispersion 
but  not  the  amplitude.  Recordings  to  the  south  show  smaller,  less  dispersive  Love  waves  as  some  paths  include 
oceanic  lithosphere  and  mantle  (see  waveforms  at  180  degrees).  Most  of  the  surface  wave  energy  is  contained  within 
the  first  pulse  and  the  amplitude  quickly  diminishes  with  successive  oscillations.  Compared  to  the  waveforms  at  180 
degrees,  those  at  350  degrees  show  an  almost  continuous  surface  wave  amplitude  of  300  to  400  seconds  duration. 
These  features  are  matched  to  first  order  by  the  current  wave  speed  model,  but  improvements,  primarily  in 
amplitude,  should  be  achievable  through  the  adjoint,  waveform  methods. 

Initial  Adjoint  Inversions 

Using  the  mismatch  between  the  observations  and  synthetic  seismograms  (e.g.,  Figure  5),  an  event  kernel  can  be 
constructed  using  the  adjoint  method  (Tromp  et  al.,  2005)  to  guide  updates  of  the  current  iteration,  wave  speed 
model.  First,  measurements  between  the  data  and  synthetics  are  made  automatically  using  the  Flexwin  tool 
(Maggi  et  al.,  2009),  which  creates  isolated  data  w  indows  and  generates  metrics  for  each  window,  that  are  then  used 
in  the  creation  of  the  adjoint  source.  Windows,  and  thus  the  adjoint  sources,  are  independent  of  any  specific  seismic 
phase  and  use  only  well  matched  data/synthetic  pairs.  Second,  adjoint  sources  along  with  the  full  synthetic 
wavefield,  are  propagated  in  reverse  time  through  the  current  iteration  wave  speed  model.  Third,  the  interactions 
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between  the  adjoint  sources  and  the  time-reversed  wavefield,  integrated  over  time,  generate  kernels  specific  to  the 
measurement  made,  e.g.,  a  simple  difference  between  the  data  and  synthetics. 
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Figure  4.  Comparison  between  tangential  data  and  3D  synthetics  for  an  event  in  the  southern  Zagros 
Mountains,  map  at  center.  Seismograms  are  plotted  against  azimuth  to  indicate  the  station 
coverage,  wave  field  coherency,  and  to  more  readily  identify  where  the  wave  speed  model  is 
deficient.  Revised  solution  is  in  black  and  the  globalcmt.org  solution  is  in  red. 
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Figure  5.  Example  of  a  multiple-event  kernel  generated  using  the  events  and  stations,  left,  for  the  Middle 
East.  Colors  for  the  shear  kernel  on  the  right  indicate  where  the  model  requires  adjustment.  A 
north  trending  slice  near  the  Persian  Gulf  demonstrates  the  kernel's  volumetric  nature. 
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Computation  of  a  single-station,  single-event  kernel,  or  a  multiple-station,  single-event  kernel,  requires  only  one 
forward  and  one  reverse  simulation.  Once  all  of  the  single-event  kernels  are  created,  multiple-event  kernels  for 
eompressional  modulus,  shear  modulus  or  density  are  produced  by  summing  the  individual  event  kernels  together. 
Finally,  it  is  these  multiple-event  kernels  that  are  used  to  update  the  wave  speed  model.  An  example  multiple-event 
kernel  for  a  test  set  of  5  events  is  displayed  in  Figure  6.  The  map  in  Figure  6,  left,  shows  the  event  and  station 
locations  used  to  construct  the  single-event  kernels,  and  the  multiple-event  kernel  for  the  shear  modulus  is  shown  on 
the  right.  A  north  trending  slice  through  the  multiple-event  kernel  near  the  Persian  Gulf  demonstrates  the  depth  and 
volumetric  nature  of  the  kernel.  Anticipated  banana-doughnut  shapes  within  the  multiple-event  kernel  in  Figure  6 
are  modified  by  a  few  important  factors:  1)  the  simultaneous  display  of  many  source-to-station  kernels,  2)  kernels 
were  computed  in  a  3D  wave  speed  model  generating  more  variable  wave  fields  and  many  off  great  circle  path 
arrivals,  and  3)  kernels  presented  here  highlight  regions  where  a  mismatch  is  found  between  the  forward  and  adjoint 
fields,  i.e.  where  the  model  requires  updating. 

Initial  tests,  as  displayed  in  Figure  6,  demonstrate  our  ability  to  generate  kernels  from  a  set  of  events  that  are  then 
required  to  update  the  eurrent  iteration  wave  speed  model.  Once  all  of  the  source  re-evaluations  arc  complete,  the 
first  in  a  series  of  model  updates  will  commence.  After  a  number  of  iterations,  the  sourees  will  be  reevaluated  within 
the  updated  wave  speed  model  to  avoid  mapping  source  errors  into  structural  features.  Iterations  will  continue  until 
the  misfit  between  the  observations  and  synthetics  has  been  sufficiently  reduced. 

CONCLUSIONS  AND  RECOMMENDATIONS 


A  seismic  waveform  database  including  source  reinterpretations  of  the  Middle  East  has  been  created.  Work  to 
characterize  events  within  and  around  the  Middle  East  has  produced  a  large  ~200  event  database  constrained  by 
teleseismic  and  regional  waveform  CMT  inversions.  All  inversion  methods,  the  CMT  inversion  and  previously  used 
regional  distance  cut  and  paste  (CAP)  (Zhao  and  Helmberger,  1994)  and  tcleseismie  (Kikuehi  and  Kanamori,  1982) 
methods,  provide  robust  estimates  of  source  depth  and  faulting  parameters.  When  used  jointly,  all  offer  an  important 
verification  of  one  another’s  result.  Event  characterization  results  from  all  inversions  agree  well  with  previously 
reported,  routine  determinations.  Our  initial  tests  using  the  3D  model  S2.9EA  have  shown  its  ability  to  accurately 
simulate  recorded  seismic  data  within  the  Middle  East  at  a  minimum  period  of  25  seeonds  and  a  dramatic 
improvement  over  a  simple  ID  wave  speed  model  at  the  same  periods.  Further,  adjoint  tomographic  infrastructure 
for  iteratively  updating  the  current  3D  w'ave  speed  model  of  the  Middle  East  is  in  place  and  test  inversions  have 
already  identified  regions  that  will  benefit  from  the  adjoint  methodology.  We  recommend,  for  the  best  possible  3D 
wave  speed  model  result  within  the  greater  Middle  East,  that  data  from  a  wide  range  of  sources  be  collected,  quality 
checked,  and  added  to  this  new  waveform  database,  particularly  seismic  data  from  remote  and/or  sparsely 
instrumented  areas. 
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